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Summary 

Plates can have more than one buckled solution for a fixed set of boundary conditions. The theory 
for the identification and the computation of multiple solutions in buckled plates is examined in this 
paper. The theory is used to predict modal interaction, which is also called change in buckle pattern 
or secondary buckling, in experiments on certain plates with multiple theoretical solutions. A set of 
coordinate functions for Galerkin’s method are defined so that the Von Karman plate equations are 
reduced to a coupled set of cubic equations in generalized coordinates that are uncoupled in the linear 
terms. An iterative procedure for solving modal interaction problems is suggested in the paper based 
on this cubic form. 

Introduction 

Modal interaction in plates is a nonlinear boundary-value problem. The plates considered here are 
loaded by edge loads in the plane of the undeformed plate. The edge loads are defined so that they 
appear in the theory as a linear boundary condition. The nonlinear part of the problem is the set of 
nonlinear partial differential equations of Von Karman plate theory. Modal interaction is concerned 
with identifying and computing all the solutions of the plate equations for a fixed set of boundary 
conditions. 

The theory for modal interaction is examined in this paper. The theory shows why problems 
with multiple postbuckling solutions require more care in computations than problems with unique 
postbuckling solutions and provides insight into a proposed algorithm for handling such problems. 

The theory for unique postbuckling solutions is often presented in the context of perturbation theory 
(ref. 1). The assumptions implied in expansions in a single perturbation parameter do not necessarily 
hold true for modal interaction problems. Continuation methods that treat the load parameter as a 
dependent variable have been implemented in finite-element codes (refs. 2 to 4). The algorithms in 
these codes assume continuity of the solutions that may lead to poor convergence properties when 
applied to modal interaction problems. 

For discrete mechanical systems, the theory (refs. f» and 6) is general enough to treat modal 
interaction problems. In the current paper, a discrete system for the plate problem is derived in 
generalized coordinates. Newton’s method is the connection between the nonlinear partial differential 
equations of plate theory and the discrete theory. The analysis in the paper starts with the linear 
form of Newton’s method. This form of Newton’s method has convergence problems for modal 
interaction problems. However, the source of the poor convergence in the linear method is clear, and a 
modification to the linear algorithm is printed in this paper, so that rapid convergence is maintained. 
The modification is an extension of previous work on postbuckling analysis using Newton’s method 
(ref. 7) for problems with isolated bifurcation points or limit points. 

Newton’s method is started by reducing the nonlinear problem to a sequence of linear boundary- 
value problems. Two difficulties arise in the application of Newton’s method. First, the linear boundary- 
value problems do not have closed-form solutions; therefore, the computation requires some kind of 
an approximate solution. The second difficulty is the hick of convergence of the sequence of linear 
problems near bifurcation points, where the different nonlinear solutions intersect. Because the linear 
boundary- value problems must be solved by approximate methods, one approach is to discretize the 
nonlinear problem from the beginning and then apply Newton’s method to the resulting set of nonlinear 
algebraic equations. 

The completely algebraic approach is not followed here because of the second difficulty of ensuring 
convergence of Newton’s method. Convergence is obtained here by converting the linear partial 
differential equations to a Sturm-Liouville problem with orthogonal eigenfunctions. The nonlinear 
problem is then discretized by a Galerkin solution with the eigenfunctions as coordinate functions. 
The resulting set of nonlinear algebraic equations with unknown generalized coordinates is uncoupled 
in the linear terms. The solution procedure can be summarized as a Galerkin solution of the nonlinear 
boundary- value problem for the plate. However, the choice of coordinate functions is not arbitrary; 
the coordinate functions are determined by the problem itself. The second, more subtle, difficulty 
of convergence is overcome directly before addressing the first difficulty, that of computing accurate 
approximate solutions for linear boundary-value problems. 



The final result in this paper is a set of nonlinear algebraic equations in generalized coordinates. 
The equations are in the form postulated by Thompson and Hunt (ref. 6) in their studies of stability 
theory. The coefficients in the algebraic equations are defined explicitly for the plate problem. The 
process of defining the coefficients suggests algorithms for numerical solutions of the plate problem. 
Newton’s method applied to the continuous formulation of Von Karman plate theory provides a direct 
connection with the theory for discrete approximations. 


Symbols 

^11^12^22^66 

Oj 

D 

T>H,Z?12,Z>22i £*66 
(E\, E 2 , £3) 

€ 3/5 Cxy 

< f,9 > 

L\2(u) 
L\2{v), L 2 2(v) 

M x , JVfy, A $xy 
Nx 1 Ny, N X y 
N\{w, w), N 2 (w y w) 
VVy, N X y , w) 
Qi 

u, v 

(u 0 ,Uo^o) 

w 

(6u, 6v , 6w) 

A 

A, 

(£ii Vi) tfii) 


stretching stiffness terms, force per unit length 

coefficient in infinite series for wq 

discriminant of cubic equation in q( 

bending stiffness terms, force times length 

residual error in plate equations for current approximation 

membrane strains 

integral of product of functions / and g over area of plate 
unit vectors 

linear terms in partial derivatives of u 

linear terms in partial derivatives of v 

linear terms in partial derivatives of w 

moment stress resultants, force 

membrane stress resultants, force per unit length 

bilinear terms in partial derivatives of w 

bilinear terms in stress resultants and partial derivatives of w 

generalized coordinate in Galerkin solution of plate equations 

in-plane components of plate displacements 

solution of linear in-plane equilibrium equation 

current approximation for a solution of the plate equations 

transverse plate displacement 

correction to current approximation of the plate equations 
load parameter that multiplies displacements on plate boundary 
ith eigenvalue of Sturm-Liouville problem 
ith eigensolution of Sturm-Liouville problem 


A subscript following a comma indicates partial differentiation with respect to the subscripted 
variable. 


Nonlinear Equations for Plate Problem 

The nonlinear theory (ref. 8) is summarized here first, and the linear form of Newton’s method is 
then applied to the three equilibrium equations written in terms of the displacement components u, 
v, and w. The linear problem is then transformed into a Sturm-Liouville problem that shows why the 
linear form of Newton’s method can fail to converge for modal interaction problems. The next section 
of the paper contains the modification of the linear form of Newton’s method using the eigenfunctions 
of the Sturm-Liouville problem. 
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The nonlinear plate equations for a specially orthotropic plate are summarized as follows: 
1. Constitutive relations: 


N x = A\\e x + A \2 * y y 

(la) 

Ny - A\ 2 e x + A 22 *'y 

(lb) 

N X y = A§§e X y 

(lc) 

Mx — -^11 ^, XX A D\2W,yy 

(2a) 

My = Di2U), xx 4* D22W,yy 

(2b) 

M X y = Dft§W' X y 

(2c) 

2. Strain-displacement relations: 


e x = u ,x A (l/2)w^ x 

(3a) 

e y = v ,y A (1/2 )wf y 

(3b) 

e xy = u,y A v x A u iX i v,y 

(3c) 

3. Equilibrium equations: 


N X)X A ^xy,y = 0 

(4a) 

N X y r X “b Ny,y = 0 

(4b) 

M. X>XX A 2A ^ r,y,xy ”b My^yy = {N x W yX -b N X yW^y) ,x A (NyW 

,y + N X yW, x ),y (4c) 

4. Boundary conditions 


The plate problems considered here have boundary arcs where u and v are prescribed. The plate 
may or may not have additional arcs that are free of membrane stresses. The boundary conditions are 
indicated schematically by t he loaded boundary C i and the unloaded boundary C 2 in figure 1. In the 
figure, the displacements on the boundary C\ are of the form 

II 

v 

e 

o' 

O 

(5a) 

v = A vi (on C \ ) 

(5b) 

where ui and vi are functions of arc length on the boundary, and the load factor A is a scalar multiplier. 
On the other hand, the stress resultants vanish on the C 2 boundary defined as follows: 

N n — N ns — 0 (on C 2 ) 

(5c) 

The boundary conditions 011 the transverse displacement w are homogeneous and correspond to simply 
supported or clamped conditions in linear-plate bending theory and are of the following form: 

w = 0 on C; iu n = 0; or w^ nn = 0 on C 

(5d) 


The boundary conditions on w could be more general, but the method of analysis is illustrated with 
less complexity by the choice made here. The load parameter (A) is also introduced in the in-plane 
boundary conditions to simplify the analysis. 

Substituting the strain-displacement relations (eqs. (3)) into the constitutive relations (eqs. (1)) and 
the resulting equations into the equilibrium equations (eqs. (4)) allows the equilibrium equations to be 
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written as nonlinear partial differential equations in terms of the displacement components u, v , and 


w. The equations, written here in operator notation for brevity, are 

L\\{u) 4 L\2{v) 4 N\(w, w) = 0 (6a) 

L 12(11) 4 L22(v) 4 N2(w, w) = 0 (6b) 

^33(^) - Nz(N x ,N y ,N xy , w) = 0 (6c) 

The linear operators, labeled Ljj, are 

Ln(u) = Ai\u, xx 4 A§§u,yy (7a) 

L\2 (/) = {A\2 + ^66)/, xy (7b) 

L22 (v) = ^22^,yt/ + AqqV }X x (7c) 

7>33(^) = Dll^^xxxx + (2Z?i2 4 4* 7^22 ^,yyyy (7d) 

The nonlinear operators are 

Wl(/,ff) = ^11 /,ite + (^12 + ^66)/,y0,xy + Awf^g^y (8a) 

^(/i #) = A22f,y9,yy 4 (^12 + ^66)/,x<7,xy 4 A^f^g^x (8b) 

Ns(N x , Ny , N xy , w) = (N x w iX 4* A/xyU^y)^ 4- (N y w, y 4 iVxyU^x) ,j/ (8c) 


Linear Form of Newton’s Method for Plate Problem 

Linear Form of Newton’s Method for Displacement Formulation 

Newton’s method starts with an approximate solution (uq, i>q, tuo) for the nonlinear system (eqs. (6)) 
plus the boundary conditions. The zeroth approximation is corrected by letting 

u = uq 4 6u , v = Vo 4- 6v, w = w$ + 6w (9) 

The linear form of Newton’s method seeks the correction (6u, 6v, 6w) by substituting equations (9) into 
the nonlinear equations (6) and dropping nonlinear terms in the corrections to arrive at the following 


linear variational equations: 

L\x( 6u) 4 L\2(6 v) 4- N\(6w,wq) 4 Ni(wq,6w) = -E\ (10a) 

L l2 (6u) 4 L,22(6v) 4 N 2 {6w, wq) 4 f> w ) = ~ E 2 (10b) 

Lm(6w) - N^(6N x ,6Ny y 6N xy ,wo) - Nz(N x o, N y o,N xy o,6w) = -E$ (10c) 

where the residual-error terms are known functions of the zeroth approximation; that is, 

E x = L n (uo) 4 Li 2 (v 0 ) 4 Ni(wo,wo) (11a) 

E 2 = Li2{uq) 4 £ 22 (^ 0 ) 4- N 2 (wq, w 0 ) (lib) 

£3 = £ 33 (^ 0 ) - Nz(N x o, N y o, N xy0 ,wo) (11c) 

N x0 = Ml [wo, I + (l/2)t«o,i + M2 [wo, I, + (l/2)to^j,j (12) 
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with similar expressions foi N y0 and N xy0 in terms of the zeroth approximation for it, v, and w. The 
corrections to the membrare stress resultants are linear in 6u, 6v, and 6w. For example, 

SNxy = A 6 Q (6u t y + 6V' X + W, x0 6w, y + 6w, x wo,y) (13) 


Linear Sturm-Liouville Theory 

The linear variational equations of Newton’s method (eqs. (10)) have variable coefficients and, in 
general, cannot be solved in closed form. However, the linear problem can be reduced to a Sturm- 
Liouville problem, which can be solved by approximate methods. The reduction to a Sturm-Liouville 
problem is achieved by defi ling the zeroth approximation for the in- plane displacements as the sum of 
two pairs of functions as follows: 


no = Xu i + u W ) (14a) 

vo = Xv L + (14b) 

The sums are defined by requiring that ui and satisfy the linear boundary- value problem, 

Lu(v-l) + L\2(vl) - 0 (15a) 

L\2( u l) + L>22{vi) — 0 (15b) 


plus nonhomogeneous boundary conditions. The functions A and Au^ that satisfy equations (15) are 
also required to satisfy the boundary conditions on u and v in equations (5). 

The functions u^o an( l <'ui0 i n th e zeroth approximation for u and v (eqs. (14)) are then defined as 


the solution of the boundai v- value problem as follows: 

£h(«ujo) + Li2{v w o) + Ni{uo,wq) = 0 (16a) 

Ll2{u w o) + L22 {v w o) + N2(uo,Wo) = 0 (16b) 

with the homogeneous bou idary conditions 

UuiO = 0 (on Ci) (17a) 

VwO = 0 (on Ci) (17b) 


Finally, the definition of the zeroth approximation for u and v as sums of solutions of two boundary- value 
problems is completed by partitioning the zeroth approximation for the membrane stress resultants 
(eq. (12) and appendix) as follows: 


N x0 = \N xL + N xw0 (18a) 

N y o = A Nyi, + N yw o (18b) 

•^YryO = A N xy L + N xyw o (18c) 

where the functions multip ied by the load parameter A and with the subscript L satisfy the boundary 
conditions on C 2 (eq. (5c)) The functions with the subscript tnO also satisfy the conditions on C 2 and 
are completely independent of the value of A. 

Splitting the zeroth approximation for u and v reduces the linear variational equations (eqs. (10)) to 
a boundary- value problem with homogeneous boundary conditions and with the load factor A appearing 
as a parameter in the following partial differential equations: 

Lji(^a) + L\2(fiv) + j\'i(6w,wo) + JVi(wg,Sw) = 0 (19a) 
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Lu{6u) + L,22(6v) + N2(6w, w o) + N2{w o, 6w) = 0 (19b) 

L33(t>u>) ~ N3(SN x ,6N y ,6N xy ,wo) - N3(N xw0 ,N yW Q,N lyw o,6w) - \N3(N X L,N yL ,N xy L,6w) = -E3 (19c) 


Eigenfunctions of Sturm-Liouville Problem 

The solution for the linear system (eqs. (19)) can be formally written as an expansion in eigenfunc- 
tions rfo, and of the Sturm-Liouville system: 

■^ll(^i) + ^12 (Vi) + 0) + N\(wg, <f>i ) = 0 (20a) 

£12(6) + L22{Vi) + N 2 (<t>i,wo) + N 2 {wo, <j>i) - 0 (20b) 

~ Ns(n xw i, Tiy W i , n X y W i, wg) — N$(N XW 0, N yw q, N X y W g, <j > l ) 
-XiNg(N xL ,N y i,N xyL ,4>i)-0 (t = 1, 2, 3, 00) (20c) 

where the functions n xw i, n yw i, and n xyw i are defined in equations (A7a) to (A7c). The notation for 
equations (20) is somewhat cumbersome, but the final result is simple. The eigenfunctions that are 
solutions of equations (20) obey the following orthogonality relation that is derived in the appendix: 

(A; — A j) J [*/*_;! A3 (Aj.£, N xy L, 0j)] dA — 0 (21) 

The orthogonality relation suggests seeking the solution of the linear variational equations as a modal 
expansion as follows: 


00 


6u = Uu,! = QiZi 

i = 1 

(22a) 

OO 

6v = v w i = ^2 qiru 
i=l 

(22b) 

00 

6w = Y,Qi ( t > i 
i=l 

(22c) 

00 

= N X w\ ~ ^ ^ Qi^xwi 
i= 1 

(23a) 

00 

6Ny = Ny W \ = Qi^lywi 

1=1 

(23b) 

00 

6N X y — N X y W i = y] Qi^xywi 
t= 1 

(23c) 

The generalized coordinates are unknown. They are determined by Galerkin’s method. The assumed 
solution satisfies equations (19a) and (19b) term by term. Equation (19c) is satisfied in the least-squares 
sense by multiplying the equation by each eigenfunction <j>j in turn, integrating over the area of the 
plate, and equating the results from each side of the equation. Because of the orthogonality condition 
(eq. (21)), the resulting equations take the form 


< —Eg <pj > 
(A -A,) 


(j = 1, 2, 3, .... 00 ) 


(24) 
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where < f,g > represents the integral over the area of the product of the functions / and g. That is, 

< f,g>= f(fg)dA (25) 

It is also assumed in equations (24) that the eigenfunctions have been normalized so that 

< <i>u N 3 (N xI , N yL , N xyL , <f>i) > = -1 (i = 1, 2, 3, . . . , oo) (26) 

Equation (24) is a canonic al form for the edge-loaded plate problem. The solution of the linear 
variational equations of N'wton’s method is completed by substituting the generalized coordinates 
Qi from equations (24) into equations (22) to determine the correction (6u,6v,6w). If the correction is 
small compared with the z *roth approximation (uq,vq,U!q), then an iterative solution of the nonlinear 
plate equations based on the linear form of Newton’s method can be expected to converge. The linear 
iterative sequence is continued by going back to equations (9), letting w\ = wq + Sw, and repeating 
the analysis through equations (24) with w\ as the zeroth approximation w 0 for w during the second 
iteration cycle. 

For the iterative procedure in the linear form of Newton’s method to converge, the corrections from 
successive iteration cycles must approach zero as the iteration continues. A quantitative measure of 
the correction 6w from any iteration cycle is the magnitudes of the generalized coordinates <7j, which 
are modal amplitudes, as determined by equations (24). The current solution for w can be expanded 
in a least-squares sense as a series in the eigenfunctions of the Sturm-Liouville problem 

00 

w = ^ a{4>i (27) 

1 

A measure of convergence is the set of ratios of the absolute values of Qi to the largest coefficient 
in absolute value in the series for w. Each q{ in turn depends on the numerators and denominators 
in equations (24). The denominators are a function of the load parameter A and the eigenvalues A;. 
Obviously, a large value for the corresponding qi is the result of dividing by a denominator that is zero 
or small in absolute value, unless the numerator is also small 
The numerators are als > linear in A and can be written as 

< £3, <f>i > = E$i - di\ (28) 

The coefficient a; in equation (28) is the same as the coefficient of fa in the expansion for w = w 0 in 
equation (27). (See appendix.) Therefore, equations (24) can be rewritten as 

Qi- (E 3 i-ai\)/(\-\i) (* = 1, 2, 3, ..., 00) (29) 

The load parameter can be prescribed to make one of the generalized coordinates, for example qj , 
vanish if the corresponding coefficient a/ is not zero, that is, 

A = E 3 j/aj (30) 

If A - A/ is also zero, <7/ can still be set to zero arbitrarily during any given iteration cycle. Then, if 
A is not too close to the remaining Aj, the iteration using the linear form of Newton’s method can be 
continued. If one of the equations is indeterminate, for example, 

q, = ( E ik ~ a k^)/(^ ~ ^k) = (0 + 0A)/0 (i / fc) (31) 

the iteration can be continued after setting equal to zero. However, the indeterminate form is an 
indication of modal interaction with other solutions to the nonlinear problem for which <7* is not zero. 
A bifurcation point or limit point, where a second nonlinear solution intersects the current solution, is 
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a limiting case of the indeterminate form in which the expansion is about an exact solution, so that 
the residual error function £3 = 0 for A = A and A = A*.. 

In cases where some of the ft determined by equations (29) are not small or are indeterminate, 
the linear form of Newton’s method must be modified either to speed convergence of the current 
approximation or to obtain convergence for solutions that intersect the current solution at bifurcation 
points. 

Modified Newton Method for Modal Interaction Problems 

The analysis of the preceding section shows that the linear form of Newton’s method for buckled 
plates may diverge or miss solutions when the load parameter A is near one or more eigenvalues A t of 
the Sturm-Liouville problem associated with the linear variational equations. A modification of the 
iteration that has better convergence properties is derived in this section. 

The linear variational equations (eqs. (10)) were derived by dropping nonlinear terms in the 
correction 6w. The exact nonlinear equations that result from substituting equations (9) into the 
nonlinear boundary- value problem (eqs. (6)) are 


Ln(Su) 4 - L a {Sv) 4 - N\(6w, wo) + Ni(wo, 6 w) = —E\ — N\(Sw,6w) (32a) 

L)2(8 u) + L 22(5v) + N2(6 w , wo) 4- N2 (wq, 6 w) = —E2 — A '2 (Aw, 6w) (32b) 

L 33 (6 w) - N 3 (6N x ,6N y ,6N xv ,w 0 ) - N 3 (N x o,Ny 0 ,N xy0 ,6w) = -E z - N 3 (6N x ,6N y ,6N X y,6w) (32c) 

The nonlinear terms in 6w are placed on the right-hand sides of equations (32) to indicate that an 
iteration sequence can be devised in which the nonlinear terms in 6w are based on some current 
approximation. The form of that iteration sequence is suggested by examining the exact nonlinear 
equations in the generalized coordinates ft. 

Nonlinear Problem in Generalized Coordinates 

In equations (14), the initial approximation for the in-plane displacements (ug^o) was partitioned 
into two sets of functions. This partitioning remains the same for the complete nonlinear problem. The 
corrections ( 6u , 6v) are further partitioned for the nonlinear problem as follows: 

Su — ^wi Uy,2 (33a) 

6v - v w \ + v w2 (33b) 

where (u^i, v w i) is the first approximation for (Su, 6v) defined by equations (22). The additional terms 

u w 2 and v W 2 in Su and 6v are defined by quadratic terms in 6w. Formally, they satisfy the differential 
equations 

Lli( u w2) + -^12(^2) + Ni(6w,6w) = 0 (34a) 

Li2{u w 2) + L22 (v W 2) + N2(6 w, 6 w) = 0 (34b) 

The nonlinear operators Ni(6w,6w) and N2(6w,6w) are quadratic in 6w. The series solution for 6w 
in equation (22c) is unchanged in the complete nonlinear formulation, and the generalized coordinates 
qi in the series for 6w remain to be determined. Equations (34) are solved in terms of the generalized 
coordinates to obtain 


OO 00 

u w 2 = 2222 (35a) 

j= 1 Jt=l 

00 00 

Vw2 = 22 22 WcVjk (35b) 

j=iifc=i 
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where the functions £ j ^ and Tfjk are defined as solutions of the following equations. 

Til (Zjk) + LuiWjk) + = 0 ( 35c ) 

Ln((jk) + L22 ( J ljk) + ^2 {4>j,4>k) — 0 ( 35d ) 

The complete corrections t > the stress resultants in terms of the generalized coordinates are obtained 
by adding quadratic terms to the linear terms already defined in equations (23). That is, 

00 00 00 

bN x = N xw \ + N xw 2 = Qi n xwi EE QjQk n xwjk (36a) 

i= 1 j=l fc=l 


00 00 00 

SNy = ^Vyu;l + Ny W 2 = ^^Qi n ywi d EE QjQk n ywjk 

i = 1 j=l *=1 


(36b) 


00 oc 00 

= N X ywl d* ^Yryu;2 — E Qi^xywi d" EE QjQkn xy wjk ( 36c ) 

i=l j=lJt=l 


The terms in the double summations are defined explicitly in the appendix. 

The final step is to substitute the partitioned form for the stress resultants into the transverse 
equilibrium equation (eq. (19c)) and to solve the resulting equation by Galerkin's method. The 
equilibrium equation is 


L^^Sw) — N${N XW \ , Ny W \ , N X y W \ , u^o) 

— Ns(N Xw0 , NywOl N X ywO) fiw) 

— \Nz(N x l, NyLi N xyL'i 
= — £3 + N%( N xw \ , ftywl » N X yw \ > ^'w) 
d" ^3 ( N xw2 * N y W 2 > ^xyw2 * ^0 ) 
d" ^3(^ r xu;2i Ny W 2i ^xyw2 t >^ w ) 

The left-hand side of equat ion (37) is identical to equation (19c). The quadratic and cubic terms in the 
generalized coordinates have been added to the right-hand side of the equation. The formal Galerkin 
solution of equation (37) is an infinite set of cubic equations in the as follows: 

(A - Aj)< 7 { - Cq 1 + C2ijkQjQk + CujtmqjqkQm (i = h 2, 3, . . . , oc) (38) 

In equations (38), repeated subscripts j , A;, and m are summed. The integrals that define the coefficients 
in equations (38) are listed in the appendix. Equations (38) are the final results of reducing the nonlinear 
plate problem to a nonlinear algebraic problem in the generalized coordinates q t . 

Equations (38) are similar to the equations in generalized coordinates postulated by Thompson 
and Hunt (ref. 6) for conservative systems. The analysis in the body of this paper and in the 
appendix gives a precise formulation for determining the coefficients in the cubic equations. The 
displacement formulation of the Von Karman equations and the choice of boundary conditions allow 
exact determination of the role of the load parameter A. The equations for the plate problem are cubic 
in the generalized coordinates and linear in the load parameter A. 

Only the real roots of the cubic equations correspond to real solutions of the plate equations. The 
number of real roots is affected by the algebraic signs of the coefficients. Therefore, in the theoretical 
analysis, much qualitative information is available when it is known how the coefficients vary with 
the load parameter A. The equations are uncoupled in the linear terms for the values of q t with the 
coefficients (A- A J. The sign of each coefficient of a linear term depends on whether A exceeds A i} which 
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in turn is a function of u;q. The coefficients and Cijk m are independent of the load parameter A 
for the edge-loaded plate. The residual-error terms Cqi are exactly the same terms that appeared in 
the linear form of Newton’s method as follows: 

C 0i = < Efrfo > =E 3i + <ii A (i = 1, 2, 3, . . . , oo) (39) 

The form of the complete set of cubic equations in generalized coordinates suggests a modification 
of the linear form of Newton’s method that will improve convergence. 

Modification of Linear Form in Generalized Coordinates 

The linear form of Newton’s method corresponds to dropping quadratic and cubic terms in 
equations (38) to obtain equations (22). This procedure breaks down when the load factor A is nearly 
equal to one of the eigenvalues (e.g., the I th eigenvalue). A modification of the linear form is to retain 
only nonlinear terms in the modal amplitude qj during a given iteration cycle. The modified form of 
equations (38) is then 


(A - A i)qi — (E 3 i 4 a* A) 4- 4 C 3i jjqqj ( i = 1, 2, 3, . . . , oo) (40a) 

The /th equation of this modified set is the cubic equation 

(A - A f)qj = (E 3I 4 ajX) 4 C^/g/gf 4 C 3 ////?/ (40b) 

If the coefficient C 3 jh is not zero, the cubic equation has at least one real solution for qj. In some 
cases, there can be three real roots. The discriminant of the cubic equation is 

D = c 3 + (36 2 /4)c 2 - (3bd/2)c - d(b 3 + d/4) (41) 


where 

c = (A - A/) /(3C 3 ///), b - C2 iii/(3C 3 jjji), and d = (E 3I - ajX)/C 3III 

When the discriminant D is positive, the cubic equation has three real roots. The real roots qj of 
the cubic equation in qj are then substituted in the remaining equations (40), which are linear in the 
remaining q t . The solution of equations (40) for any root of the cubic equation completes an iteration 
cycle in the modified form of Newton’s method. 

This solution of equations (40) can be the basis of a second iteration cycle to compute a solution 
of the complete set of nonlinear algebraic equations (eqs. (38)). The solution of equations (38) can be 
assumed as 

Qi = + s 9i (i = 1, 2, 3, . . . , oo) (42) 

where qj ' denotes the q j from a solution of equations (40) and the values of bq t are corrections to be 
determined. Substitution of equations (42) into equations (38) results in a new set of coupled cubic 
equations in the unknown <5g;. An approximate solution of this set is obtained by truncating the new 
set of equations; this truncation is accomplished by retaining only nonlinear terms in 8qj to obtain an 
updated set of equations (40); these equations are solved to complete the second modified iteration 
cycle. The modified iteration can be continued until an accurate solution or solutions are obtained 

for qi in equations (38). When the absolute values of g^ are all small, the iteration can be expected 
to converge; this convergence is expected, since the residual-error vector for any iteration cycle after 
the first is equal to the summation of the quadratic and cubic terms dropped in going from the full 
nonlinear set of equations (eqs. (38)) to the truncated set (eqs. (40)). 

This direct iterative solution of equations (38) is straightforward in theory but has disadvantages 
for actual computations. The large number of coefficients of quadratic and cubic terms are defined 
by integrals that must be evaluated numerically or in closed form. An alternative approach, which 
is equivalent to summing the quadratic and cubic terms into an updated error vector, is to update 
the zeroth approximation for wq after solving equations (40) for the first time. The current correction 
6w is computed from equation (22c). Examining the details of updating wq is beyond the scope 
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of the present paper. However, if the numerical analysis is also connected with a Galerkin solution 
for the linear eigenvalue problem, updating icq becomes a part of an iterative procedure whose rapid 
convergence makes it a practical numerical algorithm. In i his procedure, it is also necessary to compute 
only a few of the eigenvalues A 2 and the corresponding eigenfunctions fa. An equivalence transformation 
can be used in the numeric al analysis (ref. 9). The equivalence transformation is derived by going back 
to equations (22) and letting the correction 6w be a series of admissible functions gi , of which only a 
small number are eigenfunc tions fa. 

Modification for Modal Interaction 

The modified form of Newton’s method can be extended to the cases where two eigenvalues (e.g., 
Ai and A 2 ) are close together. In those cases, equations (d8) are truncated during an iteration cycle to 
contain only nonlinear terms in q\ and q*i- The pair of equations for i = 1 and i — 2 are simultaneous 
cubic equations. Real solutions for q\ and q 2 are then substituted in the remaining linearized equations 
for the rest of the q{. 

Conclusions 

Newton’s method has been applied to the nonlinear postbuckling problem for plates. The method 
reduces the nonlinear partial differential equations of plate theory to a set of simultaneous cubic 
equations in generalized coordinates. The cubic equations are uncoupled in the linear terms. The 
uncoupling is achieved by solving the linear variational equations of Newton’s method as a Sturm- 
Liouville problem. The eigenfunctions of the Sturm-Liouville problem are then used in a Galerkin 
solution of the full nonlinear plate equations to derive the set of cubic equations. 

By specifying boundary conditions on displacements, instead of on the in-plane stress resultants, 
the coefficients in the cubic equations are linear in the load parameter, which is a multiplier of the 
boundary conditions. The analysis also shows that coefficients of quadratic and cubic terms in the 
generalized coordinates of the cubic equations are independent of the load parameter. 

The special form of the cubic equations suggests a method of solution for modal interaction problems. 
The method is a modification of the linear form of Newton’s method. The solutions of the plate 
equations derived by the method can be either approximate or very accurate, depending on the number 
of generalized coordinates retained in the solution. 


NASA Langley Research Center 
Hampton, VA 23665-5225 
July 27, 1989 
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Appendix 

Definitions of Sums of Functions and Integrals That Appear in the Analysis 

Green’s Theorem and Properties of Bilinear Operators 

The analysis in the body of the paper is an examination of a formal Galerkin solution for the 
Von Karman plate equations. The integrals in the Galerkin solution make repeated use of Green’s 
theorem for integration in a plane. The form of the theorem used here is 

I '’(r* + ‘»)*A = f c * (ri + sj) • n ds — J {jz,x + & z ,y ) dA (Al) 

The area integrals are over the area of the plate indicated schematically in figure 1, and the line integral 
is over the boundary of the plate. 

The plate equations are nonlinear and contain a number of bilinear operators. The analysis in the 
paper uses the properties of bilinear operators of the general form N(f, g). The operators Ni(w, w) and 
N 2 {w, w), defined in equations (6), are bilinear operators; each is a sum of products of linear operators. 
The operator , w) is a sum of bilinear operators. The property of the bilinear operators 

that is used repeatedly here is that if a, 6, c, and d are constants, and if 

/ = afi + bf 2 and g = cg\ + dg 2 

then 

N(f,g) = acN(f\,gi) + adN{fi,g 2 ) + bcN(f 2 ,gi) + bdN(f 2 ,g 2 ) (A2) 


For example, 

N\(w, w) = Ni(wo + 6w,wq + 6w) 

Ni{w,w) = Ni(w 0 ,w 0 ) + N\(wq, 6w) + Ni(6w,w 0 ) + Ni{6w,6w) 

and 

( 00 \ 00 

i=l / t= 1 

( 00 00 \ 00 00 

i=l j = 1 / t=lj=l 

A bilinear operator on a single sum results in a quadratic form with double subscripts. The double 
sums in the quadratic forms for the bilinear operators on w = w 0 + 6w are the main factors in defining 
the solutions for u and v in the analysis of the plate equations leading to the cubic equations in the 
generalized coordinates q j. Since the operators on u and v are linear in the first two equations of 
equations (6) and the boundary conditions on u and v are also linear, the linear partial differential 
equations whose solutions determine u and v can be solved by solving equations (15), (16), (20), 
and (35) separately and superposing the results. The superposition on u and v is summarized in the 
following section. 

Definitions of Sums of Functions 

The in-plane displacements are defined in the form 
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U = \Ul + U W Q + Uu,! -t- u w2 
V = \vi + V W Q -I- + v w2 


(A3a) 

(A3b) 



The membrane strains ire summed as 


~ ^*'xL ~b CxwO “b € X u)\ “b € X w2 

(A4a) 

e y = Ac yL “b CywO “b €ywl b 6yw2 

(A4b) 

'xy = A e xyL “b G X ywQ + e xywl "b C xyw2 

(A4c) 


where 

e xL = u L y x 
CyL = VL.y 
e xyL = u L,y *b v L,x 
e xwO = u wO,x “b (l/ 2 )t^o,x 
&xywO = u wO,y “b ^u;0,x ~b U^o ? x^O,j/ 

&xwl = ^wl,x *b ^0,x^^\x 
^xywl = u wl,y “b Vwl,x “b ^0,x "b 

Individual strain terms not listed can be derived by permutations of u and v and x and y in the terms 
listed above. 

Membrane stresses follow the same notation pattern as the strains 


JV X = "b N X wO “b N X wl "b N X w2 

(A5a) 

Ny = A A lyfc + NyyjQ “b A/y-y? | ~b A/ yyj2 

(A5b) 

Nxy — ^N X yL “b N X y W § + Nyy W \ 4* Ny w 2 

(A5c) 

where 


N x L M\ e xl “b Ayi^yL'. N xwk = ^ll^xwk *b ^12^-ywk 

( k = 0, 1, or 2) 

NyL i, + ^ywk ^12^xuifc ”b ■^•22^’ywk 

(k = 0, 1, or 2) 

N X yL ^66^x7^’ 9-nd A ? X y W k = AQ§e X y%uk 

( k = 0, 1, or 2) 

Boundary Conditions 


The solution {ul,vi) satisfies the linear partial differential equations (eqs. (15)) plus the boundary 
conditions (eqs. (5)). The functions A r x / y , A^, and N X yi satisfy the stress- free boundary conditions 
on the C2 arcs. 

Since the solution (ul, >l) satisfies the nonhomogeneous boundary conditions, the solutions that 
are added must vanish on the boundary. They are defined here to vanish term by term as follows: 

u wk — 0 on arc C\ (fc = 0,1, or 2) 

(A6a) 

v wk = 0 on arc C\ (A; = 0, 1, or 2) 

(A6b) 


The stress resultants with the subscripts wk satisfy the stress- free boundary conditions on arc C2. These 
boundary conditions appear in the derivation of the orthogonality conditions for the Sturm-Liouville 
problem that is derived from the linear variational equations of Newton’s method (eqs. (10)). 
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Orthogonality Relation for the Sturm-Liouville Problem 

The orthogonality relation for the Sturm-Liouville problem (eq. (21)) follows from the boundary 
conditions, the definition of terms in the subscript wl, and Green’s theorem. In equations (22), the 
functions u w \ and v w \ are written as sums of functions that satisfy the linear partial differential 
equations (20). From the definition of the strains in equations (A4) and the stress resultants in 
equations (A5), the terms in the summation in equations (23) are also defined as follows: 


Tl xwi — + WQ x (j)i x) T Ayi (Vi,y "b w 0.,y ( t > i,y) (A7a) 

nywi = •'4l2(£i,x + ^0 ,x<t>i,x) + ^22 (Vi,y + w 0,y<t>i,y) (A7b) 

n X ywi = Aq6{Zi,v + *7i,x + w 0 ,x4>i,y + <!>i,x w G,y) (A7c) 

Equations (20a) and (20b) have the alternate form 

n xwi,x "b ^xywi,y = 0 (® 1, 2, 3, . . . , oo) (A8a) 

n xywi,x ~b n ywi,y = 0 (^ 1? 2, 3, ... , oo) (A8b) 


The orthogonality relation is derived by applying Green’s theorem to the above equations and to 
equation (20c). For equations (A8), the result is 

Ia ^ nxwi ' x + Hx V wi ^ = Jc + nx V w d) ' nC * s ~ f A ( n xwi€j,x + n xyuii£j,y) dA — 0 (A9a) 

r }j( n xywi,x + n ywi,y ) dA = ^ rij( n xywA + n ywij) ■ nc ^ a — J^( n xywi T 1j,x + n xwiVj,y ) dA = 0 (A9b) 

The boundary conditions on the solutions of equations (20) are that each eigenfunction £, and % must 
vanish on the boundary C\ and that the dot products in the line integrals in equations (A9) vanish on 
the boundary C 2 . Therefore, in the notation of equation (25) for integrals of products, 

[ ( n xwi€j,x + n xyw i£,j t y) dA = < n xw i,£j jX > + < n xyw ii£j,y > — 0 (AlOa) 

J A 

< n X ywi> Wj,x > + < Vj,y > — 0 (i = 1, 2, 3, oo; j ~ 1, 2, 3, . . . , oo) (AlOb) 

It is assumed that each of the functions fa in the solutions of equations (20) satisfies the same 
boundary conditions as w . Green’s theorem applied to any N% operator in equation (20c) then has the 
general form 

< (j)j, N$(N X i Ny , N X y, w) > = — < (N x 4 - N X yW } y) 1 4>j^x > — < (NyW^y + N X y W > 

where 

N%{N X , Ny , N X y , w) = (N X W, X + N X yW,y), X + (NyW^y 4 " N X yW, X ) y y 

Specifically, 

< fa, N${n XW i,ny W i,n X y W i,wu) > = — < (n xwl wo, x 4 n xyw iWQ, v ), <t>j yX > - < (rcyun™o,y + n X ywiWQ<x),<f>j,y > (Alla) 

< fa, Ns{N xw o, N yw o, N xyw Q,<t>i) > = - < (N xw o 0i,x 4 N xyw o fa y ),<t>j, x > - < {Nywofay 4 N xyw o<f>i iX ),<f>j } y > (Allb) 

< <t>j,N 3 (N xL , Nyl , N xy i , <t>i) > = - < (NxL<t> t,x 4 N xl ,L<t>i >y ),<t>j,x > - < {NyL<t>i,y + WryL^i,*)’ (p Jy y > (Allc) 

The orthogonality relation (eq. (21)) is proven by multiplying the jfth equation of equations (20c) 

by fa , integrating over the area of the plate, and subtracting the result from the integral of the ith 
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equation multiplied by <f)j. Considering one operator at a time in equations (20c), it can be shown first 
that 

< > ~ < > = 0 (A 12a) 

Next, by inspection of the right-hand sides of the last two of equations (All), it is apparent that 

< (f>j , A3 ( N xw q , N X yw(), 4>i ) ^ (f)j , A^3 ( NywOf ^xywO'i ) = 0 (A 12b) 

< <t>j , ^3 N xL ,NyL, N xyL , fa ) > = < fa, N3(N xL ,N yL ,N xyL . <j>j) > (A12c) 

If equations (A10) are subtracted from the right-hand side of equation (Alla) and the terms are 
rearranged, the result is 

< (j)j , A/3 {flxwif Tlywi > ywii ^ ; o) ^ = — ((j,i T ^O.x^j.x) ^ — < Wywi' (Vj,y "b ^ 

— < (£j^ ~f T)j,x ~b ^0 ,y$j,x T ^QiX$j t y) ^ 

Finally, comparison of the above equation with equations (A7) shows that the subscripts i and j can 
be interchanged with the r-sult that 

< <fij, Ns(n xw ii Tiy W , Uy W i^Ti X y W {,W{)) > — < <fii, N3(n xw j,riy W j,n X y W j, WQ) > = 0 (A12d) 

Equations (A12), based or the solutions of equations (20), are sufficient to derive the orthogonality 
relation (eq. (21)). 

Expansion of w in Terms of Eigenfunctions fa 

In equation (27), the current approximation wq for w is expressed as a series of eigenfunctions. The 
coefficients in the series are computed from the orthogonality relation (eq. (21)) and the normalizing 
equation (eq. (26)) as follows: 


^ Ns(N x L ? Nyh> A xyL ? ^0) ^ ^ <Pi i A3 


^xL' A yL 5 ^ xyL j ^ a j 4>j 

j= 1 


> = —a, 


(A13) 


The formal expansion of as a series of eigenfunctions appears again in the analysis in equations (28) 
and is repeated in equations (39). Using equations (18), the residual error £3 is written as 


^3 — £33(^11) — N^NxwQiNywQ) N xyw o,W{)) - \Ns(N x l, N y i, N xyL , wq) (A 14 ) 


Equations (A14) and (A13 are used to derive equation (28); also, the term 

^3 1 ^ 1 ^33( u, 0 ) > ~ ^ Nxw(U NywOi N X ywO> ^0) ^ 

is independent of the load i actor A. 

Coefficients of Higher Order Terms in Generalized Coordinates 

The coefficients of quadratic and cubic terms in the generalized coordinates that appear in 
equations (38) are defined as follows: 
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< <f>i, N$(N xw i, N yw \, N xyw i, 6w) > + < <fii, Ns(N XW 2, Nyw2i N X yw2-> w o) > 

00 oc 

QjQk < 4>i,Ns(n xwj i n ywj i Tlxywj » h)> 

j=z \ k=l 

00 00 

+ QjQk ^ 0t» ^3( n xujjA:i n ywjki n xywjki 0 fc ) ^ 

j=ljfc=l 

00 00 

= S C 2ijkQjQk (A15) 

j=l fc=l 


< ^3(^xtw2) Ny W 2, ^xyw2t ^0 ^ 

00 OO 00 

— ^ y y QjQkQm ^ ^3i n xwjk’> n ywjk^ n xywjk i > ^m) > 
j=l fc=l m=l 

00 00 00 

= EEE QjQkQ™-C%ijkm (A 16) 

j=l fc = l m=l 

When equations (A15) and (A16) are substituted into equations (38), the summation signs are 
suppressed to condense the notation. 


16 



References 


1. Budiansky, Bernard: Theory of Buckling and Post-Buckling ESehavior of Elastic Structures. Advances in Applied 
Mechanics , Volume 14, Chia-Shun Yih, ed., Academic Press, Inc., 1974, pp. 1-65. 

2. Riks, E.: An Incremental Approach to the Solution of Snapping and Buckling Problems. Int. J. Solids & Struct 
vol. 15, no. 7, 1979, pp. 529 551. 

3. Crisfield, M. A.: A Fast In remental/Iterative Solution Procedure That Handles “Snap- Through.” Computational 
Methods in Nonlinear Structural and Solid Mechanics , Ahmed K. Noor and Harvey G. McComb, Jr., eds., Pergamon 
Press, c.1981, pp. 55-62. 

4. Bergan, P. G.; Horrigmoe, G.; Krakeland. B.; and Soreide, T. H. Solution Techniques for Non-Linear Finite Element 

Problems. Int J . Numer. Methods Eng., vol. 12, no. 11, 1978, pp. 1677-1696. 

5. Thompson, J. M. T.; and Hunt, G. W.: A General Theory of Elastic Stability. John Wiley & Sons, c.1973. 

6. Thompson, J. M. T.; and Hunt, G. W.: Elastic Instability Phenomena. John Wiley & Sons, c.1984. 

7. Thurston, G. A.: Continuation of Newton’s Method Through Bifurcation Points. J. Appl. Mech ., vol. 36, no. 3, 

Sept. 1969, pp. 425-430. 

8. Chia, Chuen-Yuan: Nonlinear Analysis of Plates. McGraw-Hill, Inc., c.1980. 

9. Thurston, Gaylen A.; and Bains, Nancy J. C . : Solution of the Symmetric Eigenproblem AX = \BX by Delayed 
Division. NASA TP-2514, 986. 


17 





WASA Repor 

Nahonai Aeronautics and * 

Space Administration 

1. Report No. 2. Governir 

NASA TP-2943 

4. Title and Subtitle 

Modal Interaction in Post buckled Plates 
Theory 

7. Author(s) 

Gaylen A. Thurston 

9. Performing Organization Name an i Address 

NASA Langley Research (’enter 
Hampton, VA 23665-5225 


Report Documentation Page 


2. Government Accession No. 


3. Recipient's Catalog No. 

5. Report Date 

November 1989 

6. Performing Organization Code 

8. Performing Organization Report No. 

L-16573 

10. Work Unit No. ~ ~ 

505-63-01-09 

11. Contract or Grant No. 


12. Sponsoring Agency Name and Address 

National Aeronautics and Space Administration 
Washington, DC 20546-0001 


13. Type of Report and Period Covered 

Technical Paper 

14. Sponsoring Agency Code 


15. Supplementary Notes 


16. Abstract 

Plates can have more than 
for the identification and i 
this paper. The theory is i 
pattern or secondary buckl 
A set of coordinate functl 
equations are reduced to 
uncoupled in the linear te 
suggested in the paper bas 


one buckled solution for a fixed set of boundary conditions. The theory 
the computation of multiple solutions in buckled plates is examined in 
used to predict modal interaction, which is also called change in buckle 
mg, in experiments on certain plates with multiple theoretical solutions, 
ons for Galerkin’s method are defined so that the Von Karman plate 
a coupled set of cubic equations in generalized coordinates that are 
ins. An iterative procedure for solving modal interaction problems is 
<*d on this cubic form. 


17. Key Words (Suggested by Authors . s)) 

Modal interaction 
Secondary buckling 
Postbuckling of plates 


18. Distribution Statement 

Unclassified— U nlimited 


19. Security Classif. (of this report) 20. Security Classif. (of this page) 21. No. of Pages l 

Unclassified Unclassified 19 

NASA FORM 1626 oct 86 

For sale by the National Technical Information Service, Springfield, Virginia 22161-2171 


Subject Category 39 

21. No. of Pages 22. Price 

19 A03 


NAS A- Langley, 19811 





